{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Building a protein-membrane molecular system (GPCR)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Specifically we are building the $\\mu$ opiod receptor using the pdbid 4dkl from the OPM database."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "HTMD Announcement\n",
      "16-Oct-2015: New version (0.0.9) of HTMD available please update with \n",
      "\t     conda update htmd\n",
      "\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/shared/sdoerr/Software/anaconda3/lib/python3.4/site-packages/pyEMMA-2.0.1-py3.4-linux-x86_64.egg/pyemma/coordinates/util/stat.py:31: DeprecationWarning: Call to deprecated function hist. Called from pyemma.coordinates.util.stat line 31. Please use pyemma.coordinates.histogram()\n",
      "  def hist(transform, dimensions, nbins):\n"
     ]
    }
   ],
   "source": [
    "%load_ext autoreload \n",
    "%autoreload 2\n",
    "from htmd.console import *\n",
    "path = home() + '/data/building-protein-membrane/'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "---- Topologies files list: /shared/sdoerr/Work/pyHTMD/htmd/builder/charmmfiles/top/ ----\n",
      "top/top_all22star_prot.rtf\n",
      "top/top_all36_carb.rtf\n",
      "top/top_all36_lipid.rtf\n",
      "top/top_all36_prot.rtf\n",
      "top/top_water_ions.rtf\n",
      "top/top_all36_cgenff.rtf\n",
      "top/top_all36_na.rtf\n",
      "---- Parameters files list: /shared/sdoerr/Work/pyHTMD/htmd/builder/charmmfiles/par/ ----\n",
      "par/par_all22star_prot.prm\n",
      "par/par_all36_carb.prm\n",
      "par/par_all36_lipid.prm\n",
      "par/par_all36_prot.prm\n",
      "par/par_all36_cgenff.prm\n",
      "par/par_all36_na.prm\n",
      "par/par_water_ions.prm\n"
     ]
    }
   ],
   "source": [
    "charmm.listFiles()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "prj = 'mornap'\n",
    "salt = 0.15\n",
    "topos  = ['top/top_all36_prot.rtf','top/top_all36_lipid.rtf', 'top/top_water_ions.rtf']\n",
    "params = ['par/par_all36_prot.prm','par/par_all36_lipid.prm', 'par/par_water_ions.prm']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Found  segment between resid  65  and  263\n",
      "Found  segment between resid  270  and  352\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "array(['A2', 'B2', 'P1', 'P2'], dtype=object)"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "prot = Molecule(path + '4dkl.pdb')#from opm\n",
    "prot.filter('protein and noh and chain B or water within 5 of (chain B and protein)')\n",
    "prot = segmentgaps(prot,'protein','P') \n",
    "pcenter = mean(prot.get('coords','protein'),axis=0)\n",
    "unique(prot.get('segid'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "memb = Molecule(path + 'membrane.pdb')\n",
    "mcenter = mean(memb.get('coords'),axis=0)\n",
    "memb.moveBy(pcenter-mcenter)\n",
    "mol = embed(prot,memb)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "coo = mol.get('coords','lipids or protein')\n",
    "m = amin(coo,axis=0) + [0,0,-15]\n",
    "M = amax(coo,axis=0) + [0,0,15]\n",
    "mol = solvate(mol, minmax=vstack((m,M)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Bond between A: [serial 17227 resid 140 resname CYS chain B segid P1]\n",
      "             B: [serial 17837 resid 217 resname CYS chain B segid P1]\n",
      "\n"
     ]
    }
   ],
   "source": [
    "topos  = ['top/top_all36_prot.rtf','top/top_all36_lipid.rtf', 'top/top_water_ions.rtf']\n",
    "params = ['par/par_all36_prot.prm','par/par_all36_lipid.prm', 'par/par_water_ions.prm']\n",
    "mol = charmm.build(mol, topo=topos, param=params, outdir='/tmp/buil',saltconc=0.15)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "mol.view()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.4.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
